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.^ ■ Abstract 

B 

'^ [ We derive a microscopic expression for the dielectric susceptibihty x of 

O I a Coulomb glass, which corresponds to the definition used in classical elec- 

trodynamics, the derivative of the polarization with respect to the electric 
field. The fiuctuation-dissipation theorem tells us that x is a function of the 

thermal fiuctuations of the dipole moment of the system. We calculate x ^^- 
(N 
^^ ' merically for three-dimensional Coulomb glasses as a function of temperature 

a 

a 
o 
o 



and frequency. 
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I. INTRODUCTION 

The subject of this paper is the calculation of the dielectric susceptibility in Coulomb 
glasses, a term used for Anderson insulators with Coulomb interactions between the local- 
ized electrons. We consider situations deep in the insulating phase, when quantum energies 
t arising from tunneling are much smaller than the other important energies in the prob- 
lem, i.e.. Coulomb interactions and random energy fluctuations. The model also applies to 
systems in the quantum Hall regime far away from the peaks, for which the conductivity is 
exponentially small as compared to e^/h and the conduction mechanism is by variable range 
hopping between localized states. The model can be easily extended to granular metals in 
the insulating regime. 

Previous calculations of the dielectric susceptibility of the Coulomb glass have used an ex- 
pression directly obtained from the analogy between Coulomb and spin glasses [ij . However, 
the non-local character of the processes involved in Coulomb glasses makes this expression 
of the dielectric susceptibility inappropriate. Furthermore, it does not corresponds to the 
standard definition of the dielectric susceptibility. The first aim of this work is to present a 
microscopic expression of the dielectric susceptibility x valid for the Coulomb glass, and to 
apply the fluctuation-dissipation theorem to this expression. We also want to calculate x 
at very low temperatures. To this end we have to take into account that interactions, and 
specially those of long-range character, drastically change the properties of systems with 
localized states 0,^ . Most properties of these systems are affected by electron correlations, 
and such effects cannot be described by one-particle densities of states or excitations. To 
deal with complex excitations, methods were developed P-[l2| to obtain the low lying states 
and energies of Coulomb glasses. 

The paper is organized as follows: Section II introduces the Coulomb glass models used 
for the numerical calculations. Section HI presents the derivation of a microscopic expression 
for the dielectric susceptibility of Coulomb glasses. Section IV describes how the low-energy 
many-particle states are obtained numerically and gives a method to calculate the depen- 



dence of the dielectric susceptibility with the frequency. In section V, we present the results 
obtained for the dielectric susceptibility of Coulomb glasses at very low temperatures and 
the dependence on temperature, and frequency. Finally, in section VI we extract some 
conclusions. 

II. MODELS 

Our next results apply to any model of a Coulomb glass but to be definite we performed 
numerical simulations for the two most common models: the standard model with a uniform 
random potential distribution and the classical impurity band model (CIB). 

Efros and Shklovskii proposed a practical model to represent Coulomb glass problems 
with localized electronic states, which has been widely used and extended j^J^ . This model 
is represented by the standard tight-binding Hamiltonian: 

g = V^,,, + V-£^ <"--J^)<"^-^'' (1) 

^ ^ 4jre„ Tij 

where ni G {0, 1} denotes the occupation number of site i. We use rationalized units, unlike 
in most works on the Coulomb glass, because in this problem they constitute the most 
convenient choice. We will consider sites at random positions, with a density p equals to 
1, and simulate the disorder by a ramdon potential in each site 0j, uniformly distributed 
between —Wjl and W 12. Vij is the distance between sites i and j according to periodic 
boundary conditions in the sense of |]l3l in the perpendicular directions to the apphed electric 



field. Charge neutrality is achieved by a background compensation charge —K at each lattice 
site. 

The classical impurity band model (CIB) is a realistic representation of a lightly doped 



semiconductor in which the random potential arises from the minority impurities |n|. Here 
we consider an n-type, partially compensated semiconductor with donor concentration A^d, 
and acceptor concentration A^a = KN-q. The Hamiltonian is given by: 



H 



e^ I ^^ {I - ni){l - Uj) ^^{l-n 



1 >p (l-nj(l-nj) _ ^ [l-nj) \ 



Aireo \ 2 ^ nj 



where the donor occupation number Ui equals 1 for occupied donors, and for ionized donors. 
The index v runs over the acceptors and Tj,^ =| r^ — r,^ |, with r^ being the donor coordinates 
and Yy the acceptor coordinates. We chose the donor density p equal to 1 and imposed again 
periodic boundary conditions. For numerical reasons we constrained the nearest neighbor 
distance to be larger than 0.5. We used K=0.5 because the interaction effects are largest 
there. 

III. SUSCEPTIBILITY OF THE COULOMB GLASS 

We now proceed to obtain a proper microscopic expression for the dielectric susceptibility 
of Coulomb glasses. It is applicable to an arbitrary three-dimensional model of Coulomb 
glass as long as the interaction between charges goes as 1/r. 

A certain analogy between the spin glass and the Coulomb glass lead to an incorrect 
expression for the dielectric susceptibility of the Coulomb glass. In some sense the local 
spin Si in the spin glass is analogous to the site occupation rii in the Coulomb glass. If 
only rij = and rii = 1 are allowed due to strong on-site interaction, then the analogy is 
with spins Sj = 1/2. The magnetic field in the spin glass is then analogous to the chemical 
potential in the Coulomb glass. However, this analogy should not be pushed too far and 
does not apply to polarizabilities. The magnetic polarizability is the change of total spin 
induced by the magnetic field, but the electric polarizability is the change in the electric 
polarization induced by an electric field, not the change in the total occupation number 
induced by a change of a global potential as the direct analogy would have it. The basic 
difference between the two susceptibilities can be understood more clearly by realizing that 
the magnetic polarization (si)^ comes from field-induced flips of spatially fixed spins, where 
(...)t refers to thermal average. An analogous electric polarization can come from flips 
of local dipoles, but in many systems it involves a field induced displacement of charges 
(there is no magnetic equivalent to this because there are no magnetic charges). Such a 
polarizability is thus not represented by (nj)^ responding to a potential (as in [Q and |]T5[ ) 



but by ^iXi{ni)T responding to a field (as in [0]). The dielectric susceptibility x then is: 



^ = 7odE (') 

where E is the total electric field. 

Our first aim is to obtain a microscopic expression for Coulomb glasses of the classical 
definition of the dielectric susceptibility. Let us assume that we apply an electric displace- 
ment D. This will induce a polarization P equal to 

i 

where Xi is the position of site i, N is the number of sites, and (A?t,j)t is the change in 
the average occupation of site i due to the applied electric displacement. In the linear 
approximation, the change in the average occupation of site i due to a general change in the 
potential is given by: 

j -^ 

where A0j is the change in potential at site j. The partial derivative appearing in this 
expression is proportional to the local susceptibihty Xij'- 

d{ni)T eeo 



% T 



■X^J (6) 



where T is the temperature, and the Boltzmann constant ks is taken to be 1 throughout 
the paper. 

The change in potential corresponding to a uniform electric displacement is A0j = 
Dxi/eo, so the ratio between P and D is: 

gjj ~ rp]\j- Z-^^i-'X.ij^j — Xo- \l) 

To calculate the dielectric susceptibility numerically it is convenient to apply the 
fiuctuation-dissipation theorem to xo i^i order to rewrite it in terms of thermal fiuctua- 
tions of the dipole moment. Taking into account the expression for the thermal average of 



the site occupation, {ni)T, it is easy to obtain that its derivative with respect to the potential 
in j, i.e., the local susceptibility Xij^ is equal to the fluctuation in the electron occupancy of 
the two sites involved: 

Xij = {ninj)T - {ni)T{nj)T (8) 

Using this equation in expression (|^ for xo, we arrive at 

ij 
1 



:((rfV-(4l) (9) 



TN 

where d = e '^- XiUt is the dipolar moment of the sample. The dielectric susceptibility is a 
function of the thermal fluctuation of the dipolar moment. 

Our computer simulations can per force involve only systems of mesoscopic size. For the 
macroscopic susceptibility we can imagine building a macroscopic system of many mesoscopic 
cubes of linear size L arranged to flll the space. Each of these samples corresponds to a 
particular realization of the random positions and energies of the sites involved. The total 
electric fleld that a given microscopic sample feels is the sum of the applied fleld, D, and 
the induced fleld. If the applied fleld is uniform, and the polarizabilities of all the samples 
were the same, the polarization would also be uniform and the induced fleld would come 
only from the boundary of the sample. We then have (in our units): 

eoE = D-P (10) 

and get from Eqs. (P) and (|IOD: 



x = -. — -• (11) 

J- — Xo 

At very low frequencies most samples are conducting and we have an effectively uniform 

distribution of xo over the computer "samples" . At lower frequencies there is no mechanism 

mitigating the effect broad distribution of xo and the argument leading to Eq. (|TT|) fails 

and the problem becomes Clausius-Mossotti-like. A general approach to this problem for 
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random media has been given in [^,0]. Here we shall avoid the inherent complications 



of such a computation and assume that Eq. (|TT]) is approximately valid even at higher 
frequencies if we use for xo a value averaged over many computer realizations. The relation 
a = iuje means that a finite dc conductivity implies an infinite dc dielectric susceptibility. A 
proper calculation of the DC conductivity can be done by percolation in configuration space, 
but it is a difficult problem requiring huge numerical efforts so that for three-dimensional 
systems we could only consider very small samples. Our approximate calculation of the 
divergence of the susceptibility allows us to estimate the variation of the DC conductivity 
with temperature as we will see. 

IV. NUMERICAL PROCEDURE 

A. Low-energy configurations 

We calculate the dielectric susceptibility at very low temperatures making use of the 
ground state and the very low-energy configurations of the systems. With the procedure 
that we briefly discuss bellow we obtain the flrst 5.000 many-particle conflgurations and 
calculate their dipole fluctuations, Eq. (^. 

We flnd the low-energy many-particle conflgurations by means of a three-steps algorithm 
This comprises local search [^,0, thermal cycling [TTI, and construction of "neighbour- 



ing" states by local rearrangements of the charges 0,^. The efficiency of this algorithm 



is illustrated in Ref. [|I0|- In the flrst step we create an initial set S of metastable states. 
We start from states chosen at random and relax these states by a local search algorithm 
which ensures stability with respect to excitations from one up to four sites. In the second 
step this set S is improved by means of the thermal cycling method, which combines the 
Metropolis and local search algorithms. The third step completes the set S by systemati- 
cally investigating the surroundings of the states previously found. At the end we only keep 
conflgurations with a flx number of electrons, so we work with canonical ensembles. 



B. Frequency dependence 

At finite frequencies only transitions witli cliaracteristic time tu sliorter tlian tlie inverse 
of tlie frequency contribute to the susceptibility. Thus, for a given frequency, we consider 
two configurations as connected if their tjj is shorter than the inverse of the frequency, and 
we group the configurations in clusters according to these connections. 

The characteristic transition time between configurations / and J is [0], 

Tij = uJq^ exp (2 ^ Tij/a) exp {Ejj/T) /Z (12) 

In this equation, the quantity ujq is a constant of the order of the phonon frequency, Uq ~ 
10^^ s^^. The sum is the minimized sum over all hopping distances between sites which 
change their occupation in the transition I ^ J. a denotes the localization radius, Ejj = 
TLna.x{Ej, Ej) where Ej is the energy of the state J, and Z is the partition function. 

We calculate the susceptibility of each cluster through Eq. (P), assuming thermal equilib- 
rium in the cluster. The glassy nature of our systems is responsible for the existence of the 
clusters, which indicate the non-ergodicity of the systems for times shorter than the critical 
time connecting all the configurations in a single cluster. Each realization of the systems will 
be in a given cluster and will not see the other clusters. The probability to be in a cluster 
depends on the history of the system and is very difficult to estimate. In order to obtain 
averages of the susceptibility, we will assume that the weight of each cluster is proportional 
to its partial partition function, which constitutes the simplest possible assumption. The 
results are finally averaged over many different disorder realizations. 

V. RESULTS 

If we take into account all types of transitions, including the slowest ones. Coulomb glass 
behaves like a conductor and is able to screen fully as its susceptibility diverges. But small 
samples may not have excitations which carry electrons across the entire sample and produce 
nearly equipotential surfaces at the two opposite edges. So we must consider samples above 
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a certain critical size at which the steady state susceptibihty x diverges. We found that this 
critical size is about 200 sites for both models (Eqs. (0) and (|1|)). Above this size the results 
are practically independent of size. 

Fig. |1| shows the average value of xo as a function of frequency for several temperatures. 
The plots are for the standard model of size A^ = 256. The localization radius is a = 0.2, 
which is maintained throughout the paper. At low frequencies xo increases with T while at 
high frequencies it decreases with T. The reason is that at small uj hopping extends over 
many configurations, and the main effect of T is to enhance the transition rates. At large 
u) hopping is between two optimal configurations for that frequency (or even two sites) and 
the main effect of T is to equalize the occupation probabilities. This bears analogy with 
uncorrelated hopping conductivity which increases strongly with T as u; — > and behaves 
as 1/T at high frequency. 

The result for the CIB model are very similar to those for the standard model and so we 
do not show them explicitly. The results for the CIB model roughly correspond to those for 
the standard model with an effective disorder energy of approximately 3. 

As already mentioned, an accurate calculation of the frequency dependent macroscopic 
susceptibility requires the distribution f{xo), not only the average value of Xo- Fig. ^ shows 
the integrated distribution P{xo) = /o^° fix'o) ^Xo ^^i N = 256 at T = 0.01 for different 
values of the frequency, cj — > (dotted curve), u = 10'^ s~^ (dashed curve), u = 10^ s~^ 
(long dashed curve). The solid curve is a plot of the function 1 — exp{ — (xA)"} with the 
parameters A = 1.5 and a = 0.6. In the range examined this form of integrated distribution 
fits (with varying A and a) fairly well our data for all T and u. The broad character of 
the distribution indicates large mesoscopic fluctuations and shows a need to examine in the 
future the accuracy of our approximation by taking proper account of the distribution of Xo- 

We define a critical time for saturation t^ as the inverse of the frequency for which the 
value of the susceptibility is 95 % of the asymptotic value for extremely long times. We 
studied the T dependence of this critical time Tc. For all temperatures considered the values 
of Tc are extremely large, which is a sign of the glassy nature of our systems. Since Tc is 



long and close to the saturation of xo, we expect our results to be rather accurate for this 
study. Fig. |^ plots the logarithm of t^ vs. T^^/^ for four sizes of the standard model and for 
two sizes of the CIB model. The data are fitted quite well by straight lines indicating that 
a similar mechanism which gives rise to the T~^/^ conductivity ||^ is also effective in the 
dielectric susceptibility. This should of course be expected because of the close connection 
between the two properties. The square root of the slope of each straight line yields a 
characteristic temperature Tq which has been often associated with variable range hopping 
theory in a Coulomb gap. In that theory Tq is given by: 

47reo Ksa 



where f3 = 2.8 in three dimensions [0. Such a theory does not take into account corre- 



lations, so comparing our results with the theory can assess the importance of many-body 
effects. The values of f3 obtained from Fig. ^ are 0.9 ± 0.2 for the standard model and 
0.9 ±0.1 for the CIB model, i.e., about three times less than in the one-particle theory. This 
indicates the importance of correlations. The results follow the same trend we found for two 
dimensional systems by a numerical simulation study of the conductivity |^ where we also 
found P systematically smaller than was predicted by the one-particle theory, and were able 
to identify specific many-body processes. 

To evaluate, at least approximately, the frequency behavior of the macroscopic suscepti- 
bility we calculated x{^) using Eqs. @ and (|Tl]). The results are exhibited in Fig. H where 
X is plotted vs. frequency for several temperatures. The data correspond to the data for 
Xo of Fig. |l|. Note that the susceptibihty diverges in the static limit at all temperatures 
examined, but at a different rate at each temperature. Coulomb glasses will screen like a 
metal if where not for their glassy nature, which produces typical times for the divergence 
of the dielectric function much larger than measurement times. 
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VI. CONCLUSIONS 

We derived a microscopic expression for the dielectric susceptibility x applicable to hop- 
ping systems including systems where interactions are important. The expression is par- 
ticularly suitable for low frequencies. It corresponds to the expression used in classical 
electrodynamics. Si and Varma [^] have recently study the same problem in the metallic 
limit of the metal-insulator transition for two dimensional systems, and obtain than the 
static compressibility vanishes at the transition. Some previous works PIJTSI used expression 
based on an analogy between spin and Coulomb glasses. We argues that these analogies 
cannot be extended to the susceptibility. The fundamental reason is that unlike in spin 
glasses the susceptibility in the hopping systems arises from non-local processes. 

The fluctuation-dissipation theorem tells us that the dielectric susceptibility is a function 
of the thermal fluctuations of the dipole moment of the system, instead of the fluctuations 
of the charge density, result that one obtains when using analogy between Coulomb and 
Spin glasses. We calculate x numerically for three-dimensional Coulomb glass systems as a 
function of temperature and frequency. We found that x diverges as the frequency tends to 
zero. One has to consider sizes larger than a critical one of approximately A^ = 200 for the 
CIB model and for the standard model with W = ?>. 

The logarithm of the critical time for saturation varies proportionally to T~^/^, the same 
dependence as in variable range hopping. The characteristic temperature for this dependence 
is approximately equal to 0.9, a factor of three smaller than the theoretical predictions for 
the equivalent constant appearing in variable range hopping. 
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FIGURES 
FIG. 1. Averaged values of xo obtained for the standard model, plotted against frequency, 

for several values of the temperature as follows: 0.006 (•), 0.008 (■), 0.01 (♦), 0.012 (a), 0.014 

{<), 0.016 (t), 0.018 (►), and 0.02 (o).The disorder energy is IF = 2 and the localization radius 

is a = 0.2. 

FIG. 2. Accumulated distribution probability of xo foi" ^ = 256 and a; — > (dotted curve), 
Lo = 10'^ s^^ (dashed curve), cu = 10^ s^^ (long dashed curve). The solid curve corresponds to the 
fit explained in the text. Xav is the average value at each frequency. 

FIG. 3. Logarithm of critical time Tc for reaching the static susceptibility as a function of T~^' ^ 
for four sizes of the standard model, A^ = 64 (•), 128 (■), 216 (♦), 512 (a), and for two sizes of 
the GIB model, A^ = 216 {^) and 512 (t). 

FIG. 4. Dielectric susceptibility x as a function of frequency for several temperatures, T = 0.006 
(•), 0.008 (■), 0.01 (♦), 0.012 (A), 0.014 (^), 0.016 (t), 0.018 (►), and 0.02 (o).The other 
parameters and the model considered are the same as in Fig. ||. 
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